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to shear distortion in current weak lensing analysis methods. Such overall calibration 



^ ' biases cannot be detected by the standard tests such as E/B decomposition or cali- 

bration with stars. We find that the non-Gaussianity of point spread function has a 
significant effect and can lead to up to 15 per cent error on the linear amplitude of 
fluctuations erg depending on the method of analysis. This could explain some of the 
discrepancies seen in recent amplitude determinations from weak lensing. Using an 
I elliptical Laguerre expansion method we develop a re-Gaussianization method which 

' reduces the error to calibration error of order 1 per cent even for poorly resolved galax- 

ly-j I ies. We also discuss a new type of shear selection bias which results in up to roughly 

■ 8 percent underestimation of the signal. It is expected to scale with redshift, inducing 

errors in the growth factor extraction if not properly corrected for. Understanding and 
correcting for such effects is crucial if weak lensing is to become a high precision probe 
of cosmology. 
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1 INTRODUCTION 

Grav itational lensing is a powerful tool to probe the large scale matter distribution in the universe (see e.g. lBartelmann fc Schneideil 
_ ^ I2OO1I for a review). As s uch it offers the possibility of breaking the degeneracies between many cosmological parame- 
! ters (e.g. IHu fc Teemar^ Il999l) and serves as a powerful complementary method to the cosmic microwave background 
5_j ' anisotropics. Its main advantage is that it is sensitive to the dark matter directly and so avoids many of the complica- 



tions present m otner cosmological probes at low redsnitt sucn as Ka, l axy clustermg or cluster survevs. Kecent Urst detections 
of shear correlations jBacon et aljbood IVan Waerbeke et alJlioonl . iRhodes et allbnOll iHoekstra et~al] I2OO2L iBrown et all 
I2OO2L Ijarvis et alJl20n2ll have alreadv demonstrated the promise of this method. Upcoming surveys such as NOAO deep sur- 
vey {www.noao.edu/noao/noaodeep/), the CFHT legacy survey {www.cfht.hawaii.edu/Science/CFHLS/) and the Deep Lens 
Survey {dls.bell-labs.com/) will go beyond simply detecting cosmic shear but will map it over large areas of the sky. Finally, pro- 
posed projects such as DMT/LSST {www.dmtelescope.org/dark_home.html), Pan-STARRS {www.ifa.hawaii.edu/pan-starrs/) 
or SNAP {snap.lbl.gov) will map a large fraction of the sky with high depth, achieving unprecedented statistical precision. 

Advances in the statistical precision require parallel progress on the systematics errors if the information from the current 
and upcoming experiments is to be exploited fully. The effects one is looking for are extremely small, so any small errors in the 
calibration or data re duction can lead to a large error in the final result. The most commonly used methods of reduction so far 
are based on work bv lKaiser et al.l il995l) ( KSB95). Recently there have been s evera l new approa ches proposed dealing wit h 
various aspects of this problem developed bv lKaisen ll200d) . lBernstein fc Jarvil ll2002f ) (BJ02), and lRefreeier fc BaconI ll200ll) . 
What is somewhat unclear from this work is how much difference these new approaches make and what are the remaining 
errors in them. While most of these papers present at least some tests, these are often done using simplified models and do not 
attempt to quantify the error as a function of size of the galaxy or their ellipticity. The original KSB95 method was calibrated 
off of simula tions, and other authors have investigated calibra tion issues in some of the variants and extensions of KSB95 
numerically jKaiseJbooA lErben et allboOll: iBacon et allboOll) . However, the calibration of more recent methodologies such 
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as BJ02 has not yet been studied in detail, nor have the various methods been compared in the same simulation. (Many of 
the simulations in the literature are not directly comparable since they have different assumptions about the distribution of 
galaxy shapes and sizes and the PSF and noise properties.) The goal of this paper is to address some of these issues in more 
detail. We pay particular attention to the problem of point spread function (PSF) dilution in realistic PSF and galaxy shape 
models, which only affects the overall calibration of the shear (if seeing is not significantly variable across the lensing field) 
and does not show up in the systematics tests. 

Extracting shear from the images is not the only place where a bias might arise in the weak lensing analysis. Other 
possibilities include (but are not limited to) intrinsic correlations, errors in the source redshift distribution, and so-called 
selection bias, where a galaxy is preferentially selected when its shape is aligned with the PSF distortion. However, there 
is another possible bias caused by the fact that a galaxy may also be preferentially selected when aligned orthogonally to 
the intrinsic shear, because detection algorithms preferentially detect circular objects. This bias depends on the detection 
threshold and the relative angular size of the galaxy and PSF; it becomes unimportant for poorly resolved galaxies (i.e. 
galaxies much smaller than the PSF) because in this case the size and shape of the galaxy's image are determined by the PSF 
irrespective of the galaxy's physical orientation. Since it is directly correlated with the signal it is very difficult to detect it 
using the standard systematics tests such as the 45 degree rotation test (E/B decomposition) and must be corrected for using 
simulations or some other method that explicitly takes it into account. 

This paper is organized as follows. In f|^we review past approaches to PSF corrections, present the formalism for PSF 
correction based on the elliptical Laguerre expansion, and derive two approximate methods to correct for non-Gaussianity 
of the PSF. We compare various methods in the literature and estimate their biases using simulations. In JSI^s discuss a 
shear selection bias and how it can affect the amplitude of fiuctuations. Conclusions are presented in SQls-nd in the appendices 
we describe technical details on convolution in the Laguerre expansion formalism, the linear PSF method, and the KSB95 
method. 



2 PSF CORRECTIONS 

In actual measurements, we do not measure directly the moments of intrinsic galaxy images, but rather we measure the 
moments of galaxy image after convolution with the PSF of the telescope. We will denote the intrinsic galaxy image by f{x), 
the PSF by g{x), and the measured galaxy image by I{x). Our eventual objective is to take the PSF g{x) (which can be 
reconstructed by considering the images of stars in the field) and the observed image of the galaxy I{x), and construct from 
these the ellipticity e/ of the intrinsic image that can (with an appropriate calibration factor) be used as an estimator for the 
local shear. 

Many approaches to this problem appear in the literature. The common theme in most of them is to measure the moments 
of the image and PSF using some weight function: 

M^'^ = / P,{x)w{x)I{x)d^x, (1) 

where x — {x,y) is the position on the sky, w is the weight function, and the Pi are polynomials in x and y. These moments 
M^^\ and their counterparts M^^^ for the PSF, are then used to compute the ellipticity estimator. The methods mainly differ 
in their answers to the following questions: 

1. Is the weight function a circular Gaussian, an elliptical Gaussian, or some other shape? 

2. Is the weight function used for the PSF the same as that for the galaxy, or is it matched to the size of the PSF? 

3. Does the method (partially) deconvolve the galaxy image or does it directly correct the galaxy moments for the effect of 
the PSF? In the latter case, which corrections are performed using perturbation theory, which are performed by comparison 
to simulations, and which are performed using simple analytical results valid in limiting cases? 

Question^has typically been decided in favor of a Gaussian weight due to its rapid convergence to zero at large radii (hence 
the polynomial moments ar e all con vergent) and absence of singula rities, as well as gen eral mathematical convenience. The ma- 
iority of recent theoretica l lKaiae3l2D00, R. efregier fc Baconll20fllh and observational ijBacon et alJl200(l . IVan Waerbeke et alJ 
l2nonL IPhodes et aljEoOll . iHoekstra et alJi2002l ) analvses have followed KSB95 in using a circular Gaussian weight function. 
The elliptical Gaussian weight function was introduced by BJ02; although it is considerably more complicated to work with, 
a general ellipse can be a better match to an object's actual shape than a circle. This property is important for those 
methods that treat the deviation of an object or PSF from Gaussianity via perturbation theory, including our linear and 
re-Gaussianization methods. 

There is little consensus in the literature on the best answer to Question |21 The advantage of measuring the PSF with 
its own weight function, independently of the galaxy, is that the PSF weight function can be matched to the size (and, in the 
case of elliptical moments, the shape) of the PSF. This is necessary for methods (e.g. re-Gaussianization) that treat the non- 
Gaussianity of the PSF perturbatively. The disadvantage of measuring moments of the PSF with a weight function matched 
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to the PSF is that features of the PSF outside the scale radius of the weight function are suppressed. Features at radii large 
compared to the characteristic scale of the PSF are present in nearly all real-world PSFs: they are introduced as diffraction rings 
in space-based observations, and as a power law [g oc r~^^^^) tail from Kolmogorov turbulence in ground-based observations. 
If the galaxy is of size comparable to or somewhat larger than the PSF, these features can fall inside the scale radius of the 
galaxy's weight function, with potentially disastrous results for attempts to unravel the effects of the PSF on the measured 
galaxy moments (see Section f2. St . The higher-order moments of the PSF must therefore be considered in or der to take into 
account features at large radius if a PSF weight matched to the size of the PSF itself is used. Numerical te sts llHoekstra et al, 
1998 ) indicate that the KSB95 method supplemented with a pre-seeing shear polarizability calculation llLuppino fc Kaisei 



1997r works best when the PSF is measured with the same weight function as the galaxy; this is unsurprising, since this 



method does not consider any moments higher than the fourth order (see Appendix lUl. 

Finally we arrive at Question |3 how do we convert the meas ured mom ents of th e galaxy and PSF into an ellipticity 
of the intrinsic galaxy image? Most methods, including Luppino fc Kaiseil 1^97), Ho ekstra et alJ il998t) . and BJ02, have 
corrected the gala^cy moments for the effect of the PSF using either perturbation theory, simple analytical corrections, or 
some combinati on thereof. (The original KSB95 work used simulations to compute the diluting effect of the PSF on the 
galaxy.) R. efregier fc Bacoiil (|200lh developed a deconvolution method based on the circular Laguerre expansion of the PSF 
and galaxy (see Section|2i2J, writing the PSF convolution as a matrix acting on the vector of galaxy moments and cutting off 
the hierarchy of Laguerre moments at sixth order to prevent noise amplification in the deconvolution (matrix inversion) . B J02 
have investigated matrix inversion methods as well using the raising and lowering operator formalism, with particular attention 
paid to the effects of the matrix inversion on noise. In this paper, we introduce and investigate two perturbative approaches, 
expanding around the case of elliptic al Gaus sian PSFs and galaxies. The "linear" method (Appendix^ is essentially a first- 
order approximation to lRefreeier fc B acon ( 2 00lh 's matrix inverse, generalized to elliptical moments. It improves significantly 
on the same order method proposed by BJ02 in their Appendix C and should be of particular interest to those using Sloan 
Digital Sky Survey PHOTO reduction for weak lensing. The "re-Gaussianization" method (Section l2.41 treats the deviations 
of the PSF from Gaussianity perturbatively in real space, without first measuring moments. In our tests this was the most 
successful method and the only one that comes close to the precision requirements of the future surveys. The introduction of 
nonlinear calculations such as higher-order perturbation theory or direct matrix inversion can be expected to further improve 
the accuracy of the method because the full nonlinear dependence of measured moments on the PSF is taken into account, 
however one must be careful to avoid amplifying noise or creating unstable estimators. 

We organize this section as follows: adaptive moment determinati on is discussed in Section 12.11 In Section 12.21 (and in 
Appendix^ , we explore the generalization of the Laguerre expansion of lPefreeierl (|200lh and B J02 to the case of an elliptical 
weight. In Section 12.81 we explore the limitations of BJ02's analytical correction for PSF dilution; the re-Gaussianization 
method is introduced in Section [2.41 We conclude by comparing the methods in simulations ('Section 12. 5|l . 



2.1 Adaptive Moments 



A substantial literature exists on approaches using circular Gaussian-weighted moments for g and I. BJ02 introduced "adap- 
tive" moments that use an elliptical weight function whose shape is matched to that of the object. The adaptive moments of 
an image I{x) are determined using the Gaussian with the least-square deviation from the image, i.e. we minimize: 

2 



E = 



I{x) — j4exp 



— -(cc — Xo)^M ^[x — xq] 



over the sextuplet of quantities (yl,a;o,M). Defining: 
1 



w{x) = exp 

we note that minimization of E yields the moments satisfying: 



(2) 



(3) 



Xo 



J^2 xw{x)I{x)d'^x 
/h2 w{x)I{x)d2x 

g/a^ ^ xo)i{x - xo)3'w{x)I{x)d'^ 
/r,2 w{x)I{x)d'^x 



(4) 



Because of the presence of the w, the adaptive moments are weighted. [It is simpler to work with unweighted moments, but 
these have divergent noise (KSB95).] 

The ellipticity of an object is then defined in terms of the adaptive second moments of I{x): 



(M,, - Myy)/T' 



(I) 
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Table 1. The A polynomials up to fourth order. 
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T^'^ ={M^^ + Myy); (5) 
the spin-2 tensor e — (e+, ex) is called the ellipticity tensor and the scalar t'^-* is called the trace. 



2.2 Elliptical Laguerre Expansion 

This section explores aspects of the Laguerre expansion introduced by BJ02. In particular, we rework BJ02's formalism in 
a form useful for studying elliptical objects. A Laguerre expansion is constructed around a Gaussian with some symmetric 
covariance matrix M; we begin by linearly transforming the Cartesian {x,y) coordinate system to another {u,v): 

= M"^/^ ( ^ \ = — ( ^^yy + ^ \ f X 

V 

where D = MxxMyy — M'^y is the determinant of M and — D{Mxx + Myy + This provides us with the property that 

^M"^a;. The La guerre basis functions are: 

^p,{u,v) = , ^ . e-<"'+"'>/^Ap,(^.,«), (7) 
VTT-L' • p'.q'. 

where: 

min(p,q) ^ ^ 

s— 

here m = p — g is the angular momentum quantum number and the + sign is used in u ± in if m > and the — sign is 
used if m < 0. (Both signs are equivalent if m = 0.) Equations JTJ and Q can be seen to be equivalent to BJ02's equation 
(6-8) through the use of the explicit expression for the associated Laguerre polynomials. Note that Ap, is a polynomial in two 
variables with complex integer coefficients, and Apg is the complex conjugate of A^p. In Tabled the A polynomials are given 
th rough fourth ord er; in the last column we give the result in polar coordinates. We note that the A polynomials are related 
to lR,efregieil (|200lh 's polar Hermite polynomials according to Apq (r cos </>, r sin <j!>) — ffp,q(r)e™'*. 
Finally, we consider the Fourier transform: 

^^^^ ^hj (9) 
The Fourier transform of a Laguerre basis function is given by BJ02's equation (6-62); using our A polynomials, it is: 

■N 2 2 

^p,(fc) = ^ e-<''"+'''^'/"Ap,(fc„,fc„), (10) 

where = p -f g is the order of the polynomial. The variables fc„ and fc„ are given by (A;„,fc„) — M^^^a: and satisfy 
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Now any image can be expanded in the Laguerre basis functions: 

7(x) =^6p,V^p,(x;M). (11) 

p,<i 

(Note that, through an abuse of notation, we are writing tp now as though it took x as an argument.) We note that there is 
a complex-conjugate relationship among the moments bpq = 6Jp if the image / is real-valued. Furthermore, by orthogonality 
of the i/)pg basis functions, there is an inverse transform for the coefficients: 



/ i^U{x-M.)I{x)A'x. (12) 

We will find it useful to use the reduced quantities Ppq — bpq /boo- BJ02 observe that the Laguerre coefficients associated with 
the adaptively determined covariance M satisfy /3oi = /3ii = /3o2 ~ 0. One can see this directly from the least-squares residual 
principle by differentiating equation Q with respect to (j4,a;o,M) and setting the result to zero. 



2.3 PSF Corrections 

PSF correction is typically performed by correcting the moments of the measured object using the moments of the PSF. There 
are generally two cases to consider here: either the PSF moments are measured with an elliptical Gaussian weight adapted to 
the PSF shape rather than the object shape, or they are measured with the same weight function as the object. The former 
approach is used by BJ02, and we will consider it in detail here; we defer investigation of PSF moments determined with the 
same weight function as the object. 

Another issue that arises is whether to use a rounding-kernel method to symmetrize the PSF. Here we assume that the 
best-fit Gaussian to the PSF is circular; this could be achieved using a rounding kernel. An alternative is to apply a shear (to 
both the image and PSF) to make the best-fit Gaussian to the PSF circular, and then apply an inverse shear after making the 
ellipticity measurement. [The transformation of ellipticities is performed according to equation ljB17|l .] The rounding kernel 
method may be better suited to elimination of spurious power from PSF anisotropy, since a rounding kernel can be designed 
to eliminate arbitrarily many spin-2 moments of the PSF (at least in principle, see e.g. BJ02); however, rounding kernels 
increase the effective PSF size which is not desirable. Circularizing the best-fit Gaussian to the PSF with a shear does not 
increeise PSF size, but the removal of spin-2 moments of the PSF (and hence of the associated spurious power) is incomplete 
as only the bo2 moment is forced to zero by this method. 

If the galaxy and PSF are both elliptical Gaussians, the PSF correction is trivial because the covariance matrix of the 
galaxy and PSF add to yield the covariance of the observed image. In this case, we have M/ = M/ + Mg and hence: 



^ {Ml-Mq):,^-{Ml-Mq)yy 

+ {Ml -Mg)^^ + {Ml -Mg)yy 

if) ^ 2{Ml~Mq)^y 

{Ml -M3)^^ + {Ml -Mg)yy' ^ ' 

If we were using unweighted moments, this scheme would be valid for any PSF and galaxy. For weighted moments, equation 
1131 1 fails, and a correction for non-Gaussianity of the galaxy and PSF is necessary. 

BJ02's prescription (see their Appendix C) for PSF corrections is as follows: begin by applying a rounding kernel; then 
the ellipticity of the measured object e'^' is transformed into the ellipticity of the intrinsic object according to e'^' = e^'^ /R, 
where: 

Tj{l-I3g^)/{1+Pg^) 

This equation is designed to reduce to the correct form in the limit of well-resolved objects with Gaussian PSFs, and in the 
limit of any Gaussian PSF and Gaussian object (whether or not it is well-resolved, see Appendix C of BJ02); however, no 
rigorous derivation has been provided. Indeed, we can show that for non-Gaussian PSFs, the above equation may be subject 
to significant error. The measured object has intensity: 

I{x)= [ /(x')fl(a^-^OdV^ Vlg,,....„ ,^ ^" f{x), (15) 



R2 



where 



n=0 



g{x)xi^...Xi^A^x (16) 



'R2 
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is the unweighted nth moment of the PSF. (Note that formally q is infinite.) In the limit of a well-resolved object, the effect 
of the PSF is dominated by the second (n = 2) moment because the series in equation 1151 converges rapidly and the first 
moment vanishes (we assume without loss of generality that the centroid of the PSF lies at the origin). In this limit, BJ02 
showed that: 

^ = 1 ^ + OiT^/Tf). (17) 

The relation between q^x + q-yy and the Laguerre moments can be computed: 

f , E (p!g!)'''"/3p, L.^,x,Aj,,(M-i/2x)e-==^"'"''2d2a; ,, 

qij = / g{x)x,XjA^x = ^^^^^^ , ^- = M,,[l + 4/3^^^] + (18) 

J^. ' Ep,(p!9!)-mJi,.Ap,(M-i/2x)e-M-i./2d2a, P22J , V ) 

where represents both higher-order terms in I3'^2 H^nd terms from the higher-order Laguerre moments bpg' with p > 2 or 
q > 2. One would thus expect that for well-resolved galaxies, the correction for non-Gaussianity of the PSF should take the 
form: 



(l-/3«)/(l+/3(^)) 



+ 0{T^/t}), (19) 



where the numerator of the quantity in brackets is valid to first order in /322'' and to zeroeth order in the higher moments. 
We note that the B J02 correction does not exhibit this behavior; its numerator has the form 1 — 2/322^ instead of 1 -f 4/332' • 
Thus we expect that for PSFs with positive kurtosis /322' > 0, equation 11411 will overestimate the resolution factor and hence 
underestimate the shear. This expectation is verified in Section r2.5l Note additionally that the above argument does not apply 
for the less well-resolved galaxies (i.e. those for which Tg is comparable to T/); in this case, the factor of 1 -f 4/332' ^lot 
correct. The correct factor is worked out in Appendix IHl 

Several methods could be used to improve equation 1141 . One approach is to compute the correction to the galaxy 
ellipticity to linear order in the /3pq; this approach is pursued in Appendix^ Unfortunately, realistic PSFs (e.g. turbulence- or 
diffraction-limited PSFs) have significant power at large radii and hence the convergence of the Laguerre expansion is rather 
slow. This is manifested in Section [2.51 where significant errors are found from terminating the Laguerre expansion after the 
22 moment. Therefore, in the next section we will pursue an alternative approach: we fit a Gaussian to the PSF, compute 
the residuals, and attempt to correct the galaxy image for the effects of these residuals. Then we apply equation 11411 to this 
corrected galaxy image. As shown in the simulations (Section 12.511 . this method yields good performance for non-Gaussian 
PSFs and galaxies, except for the most elongated galaxies. 

2.4 Re-Gaussianization 

Let us suppose that the PSF g can be approximated by a Gaussian G with some covariance Mg: 

We may then define the residual function e(x) = g{x) — G{x). It follows that the measured image intensity will satisfy 
I = g®f = G®f-\-e®f, where ® represents convolution. This can be rearranged to yield: 

G(g)f = I~e(g)f. (21) 

This equation thus allows us to determine the Gaussian-convolved intrinsic galaxy image /', if we know /. At first glance 
this does not appear helpful, since if we knew / it would be trivial to determine G ® f . However, / appears in this equation 
multiplied by (technically, convolved with) a small correction e, so equation 12H may be reasonably accurate even if we use 
an approximate form for /. The simplest approach is to approximate / as a Gaussian with covariance: 



/o -i exp(-lx^[Mf']-^. 



27r^detM^"' 

M!-p ^Mz-Mj, (22) 

where M/ and Mg are the adaptive covariances of the measured object and PSF, respectively. Then we can define: 

I' = I-e®h{^G®f). (23) 

The adaptive moments of /' can then be computed, and the PSF correction of equation 1141 can then be applied to recover 
the intrinsic ellipticity e'-'^-'. 
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Figure 1. The error in reconstructing the shear for the case of a circular Gaussian PSF and elhptical Gaussian galaxy. We have plotted 
the shear calibration error S'y/'y. Note that the adaptive moments methods are exact for this case. The KSB95 method is not exact but 
gives good performance (calibration errors of order 1 per cent for most galaxies) even when the galaxy is not well resolved. 



0.1 



0.05 



o 
u] 



O -0.05 



-0.1 



(a) a2/of=0.2 



1 


1 1 1 1 1 

NC 


1 




BJ, G — 






LIN 






KSB 






1 





1 



Shear Calibration Errors fpryarious PSF Corrections 
(b)o2/o2=1.0 

0.1 



o 

UJ 



CO 



0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 



0.05 - 



- O -0.05 - 



-0.1 



1 1 1 


1 1 1 

NC 


1 




BJ, G — 






LIN 






KSB 




















o 

lii 



Intrinsic Ellipticity, e' 



(f) 



0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 

.(f) 



(c) 02/02,2.0 



0.1 



0.05 - 



O -0.05 - 



-0.1 



1 1 


1 1 1 1 

NC - 


1 




BJ, G - 






LIN 






KSB 














0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 



Intrinsic Ellipticity, e* ' Intrinsic Ellipticity, 

Figure 2. The error in reconstructing the shear for the case of a circular Gaussian PSF and elliptical exponential galaxy. We have plotted 
the shear calibration error S-yf-y. Note that for the well-resolved galaxy (a), both BJ02 and linear methods provide good performance; the 
error in the linear method grows substantially for the poorly resolved galaxies in (b) and (c). The no-correction method underestimates 
ellipticities in all cases, with the situation worsening for the poorly resolved galaxies. 



There are several "knobs" that can be adjusted on the re-Gaussianization method. First, there is the choice of Gaussian 
approximation to the PSF, equation I2UII . We have chosen to set Mq to be the adaptive covariance of the PSF, i.e. the 
covariance of the best-fit Gaussian. However, rather than using the normalization of the best-fit Gaussian, we have normalized 
equation 1201 1 to integrate to unity. While this increases the overall power J (e^) of the residual function, it yields J e = 0, 
which ensures that for well-resolved objects (i.e. objects for which the PSF is essentially a 5-function), the "correction" e(8)/o 
applied by equation H23^ does not corrupt the image /. One could also try a more sophisticated function fo than a Gaussian, 
perhaps based on a partial deconvolution of the PSF; we have not explored the possibilities here, but it is clear that an 
attempt to approximate / will have to make some assumptions regarding the behavior of / in the regions of Fourier space 
where the PSF has little power (i.e. where the optical transfer function is close to zero). 



2.5 Simulation 

We may test the PSF correction algorithms via simulation for common categories of galaxy and PSF. The simulation procedure 
is as follows: 

• Construct an intrinsic image / of the test galaxy; 

• Construct a sample PSF g; 

• Convolve the test galaxy and PSF, then pixelize this to yield a "measured object" 7; 

• Compute the adaptive moments of the test galaxy /, and measured galaxy /; 
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Figure 3. The error in reconstructing the shear for the case of a circular PSF generated by atmospheric turbulence and an elliptical 
Gaussian galaxy. We have plotted the shear calibration error S^/'y. For well-resolved (a) and moderately-resolved (b) galaxies, the linear 
method substantially reduces errors when compared to the no-correction method. For poorly resolved galaxies, the linear method over- 
corrects for PSF non-Gaussianity. The PSF non-Gaussianity correction of BJ02 has the incorrect sign in all three cases. Note that the 
scale is different from that of Fig. I^J 



• Apply the PSF correction scheme to the moments of / to yield the "recovered" ellipticity of the test galaxy; 

• Compare the true and recovered ellipticities e'-'^-' of the test galaxy; 

• Propagate the ellipticity error through to an error in the shear estimate. 

Our simulations here will compare five methods of recovering the initial ellipticity. The re-Gaussianization ("G") method 
is that described in Section 12.41 The linear ("LIN") method follows our prescription at the end of Appendix^ The "BJ" 
method uses BJ02's resolution factor from Appendix C, equation l|14|l . instead of ours, equation ijBlBfl . [BJ02 discuss several 
methods of correcting for PSF dilution, including a matrix inversion method which from theoretical expectations should do 
better than the resolution factor approach. The resolution factor used here is the same one that is used in the cosmic shear 
study of Jarvis ct al. (2002) and can be implemented on Sloan Digital Sky Survey PHOTO reductions.] The no-correction 
("NC") method computes ellipticities based on equation (1131 and hence is correct only for Gaussian galaxies and PSFs. 
Included for comparison is the "KSB" method, which uses the KSB95 method with pre-seeing shear polarizability estimated 
by the method of LuDoino & Kaiser (1997). Note that t here are several versions of the KSB method in the literature, a 
comparison of some of these is provided bv lErben et alJ j200lh . (For details of our implementation of KSB95, see Appendix 

10) 

As our objective is to correct for the dilution and smearing of th e galaxy image by the PSF, we do not include noise in 
these simulations; thus noise-rectification biases such as iKaiseJ i2000l) 's "noise bias" (termed "centroid bias" in BJ02) are not 
addressed here and a separate correction for these must be made. The bias in the PSF-corrected ellipticity e introduced by 
noise is generally given by the Taylor expansion; 

where li is the fiux in pixel i. Here the leading order term scales with the significance f as l/v'^, and the higher-order 
terms scale as l/v^, etc. For galaxies detected at high significance, the dominant noise bias contribution is from the 

leading l/u^ term (shape measurement of low-significance galaxies is unstable an d not recomm ended). Present approaches to 
the evaluation of the noise bias include analytic solutions for Gaussian objects jKaiseJboOol) and analytic scaling relations 
whose coefficients can be determined through simulations (BJ02). There is also a Monte Carlo approach, valid for any object 
profile and any PSF method: if we add artificial Gaussian noise to the measured image with covariance (^{AliAIj} (where 
the parameter <^ 1), and then take an ensemble-averaged change in measured ellipticity over many realizations of the 
artificial noise, we find the ensemble-averaged change in measured ellipticity E = ^(Ae). This allows us to determine (Ae), 
and thus correct the measured ellipticity for noise bias: e i-^ e — (^^^E. Note that the statistical variance in e resulting from 
the Monte Carlo method after A'^ realizations of the artificial noise have been computed is (Ae^)/A whereas the statistical 
noise in the galaxy shape measurement is (Ae^); hence the Monte Carlo convergence rapidly becomes a negligible contributor 
to the overall statistical uncertainty. Thus we can determine and remove the noise bias effect on any PSF correction method 
by evaluating Eq. (1241 . provided only that the significance v is large enough to treat the noise as a perturbation. Since the 
0{l/i'^) noise bias is straightforward to remove and in any case the details of its removal process do not depend strongly on 
the PSF correction scheme, we have presented the simulations here without noise to illustrate the purely PSF-induced effects. 
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Table 2. PSFs used in the simulations. The adaptive (i.e. iterated until /3ii = 0) radial 27ith moments /3nn characterize the non- 
Gaussianity of the PSF. Note the slow convergence of the Laguerre expansion for the diffraction-limited PSF. 



PSF type 


equation 


adaptive width, a 


/322 


/333 


/344 


/355 


Gaussian 


eq. 




u 


0.000 


0.000 


0.000 


0.000 


turbulence-limited 


eq. 


rm 


0.45epwHM 


+0.046 


+0.012 


+0.010 


+0.006 


Airy (diffraction-limited) 


eq. 




0.41A/D 


-0.043 


+0.023 


+0.039 


+0.022 


quartic-Gaussian 


eq. 


03 


cr 


+0.046 


0.000 


0.000 


0.000 



Since the end result of the shape measurements is to produce a shear estimate, we show the calibration error 57/7 in the 
shear estimator as a function of object and PSF parameters. For the adaptive moment methods (no correction, linear, BJ02) 
that yield ellipticity estimates for individual objects, we proceed as follows: consider the (unweighted) shear estimator: 

where e is the measured ellipticity. For objects of a fixed ellipticity this has calibration: 

97- 
97- 



^ ^ (26) 



Now if we use the transformation law for ellipticities under infinitesimal shear, i--> + 27^ — 2eae/37/3 [see our equation 
(41) in the 5 — > limit, and note that the infinitesimal shear 7 differs from that of the distortion shear 5 by a factor of two], 
we can azimuthally average to yield: 



974 
974 



7=0 



This is the equation we use to determine the calibration for the no correction, linear, and BJ02 methods. In the plots (Figs, 
w'e show the error ^7/7 = g^|7=o — 1- 

For the KSB95 method, equation 12711 is not applicable because the method does not produce an estimator for the adaptive 
ellipticity. Rather it produces a non-adaptive ellipticity and a pre-seeing shear polarizability tensor P^. The shear is then 
estimated via: 

1. = [{P'r'U^ie^}. (28) 

For an ensemble of galaxies of the same ellipticity but rotated by different angles, (P^) becomes isotropized and is equal to 
TrTr P^. We can find (e^) by an argument similar to that leading to equation H27|l : the result is: 



974 
974 



Tr P^ V ^ ' de e 

7=0 



This is the calibration equation that we use for the KSB95 method. 

We use several different forms for the PSF in our simulations; these are given in Table |5| The first is the Gaussian form, 
with PSF g{r) and optical transfer function (OTF) g{k) given by: 



5(fc) = 2ng(k) = / e'"' '' g{x)d^x = e'" . (30) 



g{x)a X — e '. 

(We define the OTF g{k) to be 27r times the Fourier transform g{k) of the PSF so that g{Q) — 1.) The Gaussian, while simple 
to analyze, is not a good approximation to the PSF for real experiments. For ground-based experi ments, the re solution is 
generally limited by atmospheric turbulence in the inertial range; the OTF for this case is given by (see lKaiseil200nl Appendix 
A): 

^ISB-)'^^(^ii)^ 

-g{k) =exp[-(fc/fco)'/'] , (31) 

where fco = 2.92/Sfwhm and the seeing Spwhm is the full-width half-maximum width of the PSF; typically this is of order 
1 arcsec. [The series for g{r) is convergent for all radii, nevertheless the summation is numerically unstable at k^r S> 1 and 
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Figure 4. The error in reconstructing the shear for the case of a circular quartic-Gaussian PSF [see equation I3iil ] and an eUiptical 
Gaussian galaxy. We have plotted the shear calibration error Sj/'y. For well-resolved (a) and moderately-resolved (b) galaxies, the 
linear method substantially reduces errors when compared to the no-correction method. For poorly resolved galaxies, the linear method 
over-corrects for PSF non-Gaussianity. The PSF non-Gaussianity correction of BJ02 has the incorrect sign in all three cases. 
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Figure 5. The error in reconstructing the shear for the case of a circular turbulence-limited PSF and an elliptical exponential galaxy. 
We have plotted the shear calibration error (S7/7. Note the significant loss of accuracy of all the methods except for re-Gaussianization 
(G) for poorly-resolved galaxies. 



for computational purposes it is better to take the Fourier transform of g{k).] For space-based experiments, the telescope is 
usually diffraction-limited and hence the OTF is the autocorrelation function of the aperture. If the aperture is circular with 
diameter D, this is: 



Ji{kor/2) 



^ arccos(fc/fco) — (fc/fco)-\/l — k^/k'g , k < ko 
, k>ko 



(32) 



where ko = 2nD/X and A is the wavelength of observation. We do not pursue an analysis of the effects of chromaticity, but 
note that (particularly for the diffraction-limited PSF, which has a steeper wavelength dependence than the turbulence- limited 
PSF) the effect may be important. 

Simulation results are shown in the figures for three cases: a well- resolved galaxy {ag/a'j — 0.2; we remind the reader 
that subscript g is of PSF and / for galaxy), an intermediate case {a^/a^j = 1.0), and a poorly resolved galaxy (ag/a^ — 2.0). 
[Here = VdetM is proportional to the area covered by the object.] The simple case of a Gaussian galaxy and Gaussian 
PSF is shown in Fig. The adaptive moments methods are exact for this case. The KSB95 method, while not exact, has 
relatively small calibration error except for the very elongated galaxies. 
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Figure 6. The error in reconstructing the shear for the case of a circular turbulence-Hmited PSF and an elHptical de Vaucouleurs 
profile galaxy. We have plotted the shear calibration error (S7/7. Note the significant loss of accuracy of all the methods except for 
re-Gaussianization (G) and KSB for poorly-resolved galaxies. 
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Figure 7. The error in reconstructing the shear for the case of a circular diffraction-limited PSF and an elliptical exponential galaxy. 
Severe loss of accuracy occurs for most of the methods if the galaxy is poorly resolved (c); fortunately, diffraction-limited PSFs occur 
only for space-based surveys in which the PSF size is usually small compared to the size of the galaxy, hence the comparatively large 
errors exhibited in (c) are unlikely to be a problem in practice. 



In Fig. |5| we test the ability of these methods to handle non-Gaussian galaxies by showing simulation results with a 
circular Gaussian PSF and an elliptical galaxy with exponential (oc e~^^^°) radial profile. The exponential profile has a 
positive kurtosis /322 ~ 0.17. Note that both the BJ02 and linear methods improve upon the no-correction method; BJ02's 
performance here is superior to that of the linear method, presumably because BJ02 works to all orders in f3'^2 ■ Once again, 
the KSB95 method is accurate at the 1 per cent level except for very eccentric galaxies. 

In Fig. 121 we display a simulation using an elliptical Gaussian galaxy and a non-Gaussian PSF. The PSF is the Fourier 
transform of an e"'-^*^^''^ function, as appropriate for observations through a turbulent atmosphere with Kolmogorov scaling 
jKaiseJboOol) . This PSF has a kurtosis of p':f2 ~ 0.046; it is important to note, however, that the PSF drops to zero slowly 
at large radii and hence it possesses higher moments, e.g. /Sgl' ~ 0.012 (see Table|5J that are not taken into account by either 
the BJ02 or linear methods. We can demonstrate that the error from the linear method is due to these higher moments by 
repeating the simulation using a PSF of the "quartic-Gaussian" form: 

.(0 ^(1 + 0.046. ->^-y/-^ + ^ )e--/---, (33) 

which has the same kurtosis as the turbulent PSF but no higher-order moments. The results using this form for the PSF are 
shown in Fig.|l] Finally, we have shown simulations where both the galaxy and PSF are non-Gaussian: one with a turbulent 
PSF and exponential galaxy (Fig.|^, one with a turbulent PSF and de Vaucouleurs profile (Fig.jSJ and one with a diffraction- 
limited (Airy) PSF and exponential galaxy (Fig. |7J . Of the methods described, only the re-Gaussianization method retains 
high accuracy in these cases. This is primarily because of the sixth-order and higher moments of these PSFs, which are ignored 
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by the BJ02 and fourth-order linear methods. For well-resolved objects, the KSB95 method avoids these problems because its 
non-adaptive weight function allows it to "see" the turbulence-induced tails and/or diffraction rings in the PSF. For the less 
well-resolved objects, KSB suffers from loss of higher-moment information, although it is apparent from the simulations that 
the problem is not as serious as for the BJ02 or linear methods. 

The simulations indicate significant (at least several percent) calibration errors for all of the currently used methods 
when the PSF is non-Gaussian and has size comparable to that of the galaxy. Both the KSB95 and BJ02 methods yield an 
underestimate of the shear for turbulence-limited PSFs, with the error for poorly resolved objects being of the order of 5 per 
cent for KSB95 and 10-15 per cent for BJ02. The errors have the opposite sign for a diffraction-limited PSF, however the 
small size (typically of order 0.1 arcsec) of diffraction- limited PSFs helps to reduce the error in this case (compare Fig.EK"C)- 
The re-Gaussianization method introduced here shows promise in reducing these calibration errors to the 1-2 per cent level 
even for poorly res olved galaxies. It has the best performance of all the methods studied here. While we have not simulated 
iRefregier fc Bacoiil 1^01) 's matri x inversion method , we n ote for comparison that in realistic Monte Carlo simulations (i.e. 
including noise and pixelization'). IRefregier fc BacorJ (l200lh found calibration errors of —3 ± 4 per cent in e+ and ± 4 per 
cent in ex using a turbulent PSF with cTg/a^ = 0.32. This would thus be comparable to our tests of well resolved galaxies, 
for which KSB95 and re-Gaussianization methods are accurate to within 2 per cent. It would be interesting to investigate 
the behavior of the matrix inversion method at lower resolution where the other methods except for re-Gaussianization show 
more significant error. 



3 SHEAR SELECTION BIAS 

Even if the PSF correction is perfect, an estimate of the shear can be biased if the obje ct detection a lgorithm preferentially 
selects galaxies elongated in one direction over galaxies elongated in another direction. iKaiseJ i200Cil and BJ02 investigate 
such a selection bias that results from asymmetries of the PSF and offer methods for correcting it. However, there is another 
selection bias that exists: suppose that we are looking at a part of the sky where the gravitational shear elongates objects 
in the x-direction (i.e. 7+ > 0). Then background galaxies intrinsically (i.e. pre-shear) elongated in the x-direction appear 
more elliptical than otherwise identical background galaxies elongated in the j/-direction (although both have the same area, 
flux, and surface brightness). Many object detection algorithms, such as significance cutoffs, will preferentially detect the 
circular objects in this situation, thus the catalog of background galaxies that results will preferentially contain galaxies 
aligned orthogonal to the gravitational shear. This results in underestimation of the lensing si gnal. We will call this effect 
"shear selection bias" to distinguish it from the "PSF selection bias" discussed bv lKaiseil (|200fi l and BJ02. 

Unlike PSF-related biases, shear selection bias is perfectly (anti-)correlated with the shear signal. This means that shear 
selection bias takes the form of a calibration error, with the fractional error depending on the object detection algorithm and 
on the underlying population of galaxies. As previously mentioned, calibration errors are difficult to detect using, for example, 
E/B decomposition of the measured shear field, which is insensitive to a position-independent calibration error. 

We organize this section as follows: in Section [3. II we develop a formalism for relating the galaxy selection probability to 
the selection bias. In Sections l3.2l and l3.3l we apply this formalism to the problem of evaluating shear selection bias for galaxy 
detection algorithms based on a statistical significance threshold. We show in Section 13.41 that for this class of algorithms, 
shear selection bias is an effect of order several per cent. Suggestions for overcoming the shear selection bias are offered in 
Section 13.51 



3.1 Effect on Measured Shear 

Here we analyze the effect of shear selection bias on measurements of gravitational shear. Following KSB95, we will use 
lowercase Roman letters for vector indices and Greek letters for spin-2 tensor (i.e. traceless symmetric matrix) indices. To 
avoid confusion, we call an object with two tensor indices a "hypertensor" . 
A shear estimator 7 will generally take the form: 

(7+,7X) = ^^T-^^»44,e^), (34) 

where the sum is over the background galaxies in the region of sky where the shear is to be estimated, ej is a measure of the 
shape of a galaxy, wj are the weights for each galaxy, and T is a calibration factor. In order to have the proper first-order 
behavior, we must have: 

rs%^r^('-W), (35) 



w)d \ djl^ 

where the d subscript indicates that the average is a weighted average taken over the detected background galaxies. Because 
shears observed in weak lensing experiments have 7^1, there is typically no need for correct 0{'y^) and higher behavior. 
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In the BJ02 formalism, the shape ej = e.j, where ej is the galaxy's ellipticity. In the KSB95 formalism, ej is a shear 
estimator computed according to: 



FT-ism— li/3 5\ 
\P ] SP ) 



where is the circular-weighted ellipticity: 

(e+,eX)= ^ 



OO /"OO 



(36) 



(37) 



oc >^ — oc 



p is the unweighted quadrupole moment of the PSF, and the hypertensors and are the smear and shear polarizabilities, 
respectively. The weight function W{r^) is usually taken as a Gaussian. The smear polarizability is chosen to eliminate the 
PSF bias (i.e. the change in measured galaxy ellipticity due to elongation of the PSF), and the shear polarizability is chosen 
so as to achieve the correct response, i.e. equation H35|l with T = 1: 



[P'' 



(38) 



Methods fo r computing P"^ can be found in terms of Laguerre coefficients in our Appe ndi x IHl and in terms of real-space 
variables in lHoekstra et ahl il998l) 's Appendix A [these methods are based on K SB95 andlLupp ino fc Kaiserl il997l) . however 
KSB95's shear polarizability contains algebra errors that have been corrected bv lHoekstra et ctlTiligg^ ]. 

The effect of shear selection bias can now be analyzed. If we have an ensemble of background galaxies gravitationally 
sheared by a small amount 7, then the shear estimated according to equation 1341 is: 



1 

T{wV) 



_ 1 (e^wP) 
^ f {wV) 



7=0 



+ 



(39) 



where V indicates the probability of a galaxy being detected. (Detection is a stochastic process due to photon counting noise. 
Assuming PSF anisotropics and camera shear are correctly removed, the terms proportional to e^P vanish. This yields: 



T{wV) [\a7'3 



wV 



1 

T 



d£_ 



\ ^ / -D ^™ \ 



13 / 7 
wd 



(40) 



where we have defined the transfer hypertensor ^. Note that if the estimator's response is correctly calibrated in the absence 
of selection biases [i.e. equation (I35II is satisfied], then we have: 



and 

{7°) = 7° + ^ (2"/3)„d7'', 



(41) 



(42) 



in which case H represents the calibration error resulting from shear selection bias. Thus if we can determine H, or even its 
average value, the shear selection bias is overcome. Below, we examine the transfer hypertensor for several object detection 
algorithms, and find that H is strongly dependent on the details of this algorithm. 



3.2 Significance Threshold Algorithms 



The simplest object detection algorithm is the significance threshold algorithm. This method convolves the measured image 
intensity I{x) with some detection kernel Ki to yield an convolved intensity J: 



J{x) 



/ I{x')Ki{x ~ x')d^x' 
Jr2 



This quantity has an expectation value: 



J{x) = {J{x)) 



R2 



Ki{x — x")g{x" — x')<i^x" 



R2 



A^x' = / f{x')K{x — x')d^x' 
Jr2 



(43) 



(44) 
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where / is the intrinsic intensity of the galaxy, g is the PSF, and K is the convolution of Ki and g. The two-point correlation 
function of the noise is simply the autocorrelation of the kernel: 

= {5J{x)SJ{x + v))=n I {z)Ki {y + z)d'z, (45) 



where n is the noise variance per unit area. We will require Ki to be normalized so that ^(0) — 1. Then an object is "detected" 
if at least one point within it has J{x) > vo, where is the cutoff significance. Note that with our normalization, we can 
think of the point x as being J(x) standard deviations above background. 

Our principal objective here is to determine the probability V of the object being detected. The definition of "detected", 
and hence the determination of V, is a delicate issue. If we simply assume some intrinsic intensity f{x) peaked at a; = 
and ask whether J(0) > i^o, we will underestimate V because J is a random field and it is unlikely that will be an exact 
maximum of J. However, the use oi v = J(0) as a detection statistic will serve as a first approximation to the shear selection 
bias, and due to the complexity of the problem it is unclear whether substantial gains can be made through more detailed 
analytical investigation. 

The probability of a galaxy being detected is given by the cumulative Gaussian distribution function: 
T = nv>v,) = ]^+'^-eri'^-^, (46) 
where v = J(0) is the expected significance of the galaxy. Defining the function: 



^, , d\nV V2^ 
1 (u; vo 



^ exp 



(47) 



dlnv i + ierfi^ 

we can then write the transfer hypertensor as: 

c d adlnV dlnu ^ d\nv 

We can compute the partial derivative dlnu/dj'^ using integration by parts: 
dlnu _ ^J^,f{S,x)K{~x)d' 



^., k2^A^fi^)]K{^^)^'x ^ ^^ J^,f{x)x.^K{~x)A^x 
' J^,fix)K{-x)d'^x ^ J^,f{x)K{-x)d^x ' 



7=0 

where the mixed hypertensor S^^ has components: 

and we have used the tracelessness of to eliminate the term. In the special case where ii" is a Gaussian, K{x) oc 
exp(— i|a;p/(7^), this reduces to: 

91n^_ 1 /j,,<5;^x,^./(x)e-l-lV2.irf2^ ^ iJ^,\x\^f{x)e-^-^'^''^ld'x_ 

-6/3. (51) 



dl^ /j^, /(x)e-l"l'/2-i-d23. aj, /^^ /(x)e-l"l'/2-id2a; 



For a single detection kernel Ki, equation I51II describes the transfer hypertensor. However, the single-scale detection algo- 
rithms considered here are not optimal for detecting objects substantially larger or smaller than the size of the detection 
kernel. For this reason, most authors have used multiple detection kernels, which we consider next. 



3.3 Multiscale Detection Algorithms 

The commonly used KSB95 detection algorithm (see, e.g. KSB95. [Hoekstra et ahlFiggsIl detects objects by running a signifi- 
cance threshold algorithm on several smoothing scales, i.e. with several functions Ki. We will analyze here only the case of a 
Gaussian PSF of standard deviation ag and a Gaussian smoothing filter of adjustable standard deviation ad- Then we have: 



1 -lr.|2/2„ 



Ki{x)= 5=e 

K{x)= '^d (52) 
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T-factors for Shear Selection Bias 




4 6 
Significance, v 

Figure 8. The T-factor T(P; !^o) = rflnP/dlnP. This factor represents the (logarithmic) increase in probability of detection V if we 
increase the significance of the galaxy. The shear selection bias increases in proportion to an appropriate average of T [see equation 1481 ]. 
Note that T — > as P — > oo; this means that there is no shear selection bias for objects far above threshold, as expected. 



where a% = a^ + a^. An object is detected if it is detected with any of these filters. We must therefore evaluate equation l|51|l 
at the most significant scale, i.e. where z> is maximized with respect to cr^. We determine z/ by integration: 



u = J(0) = 



The maximum siginficance occurs where dv jdod ~ 0, that is where: 



l^e-l-l'/^-^d^a 



/^,/(x)e-l-l^/2'^Kd2x 



2 



2 



(53) 



(54) 



(Remember that ok depends implicitly on ad-) In the special case of a circular Gaussian object / oc exp(— i|a;p/(7^), the 
left-hand side of equation H54^ becomes 2/{aJ^ + oT'^)', this gives a solution to equation H54|l of: 



2 



1-^ 

2 



(0-9 + o-d) 1 1 



2 



(Jg+af. (55) 

(We ignore the imaginary solution — —cr^-) This means that if the galaxy is a circular Gaussian, the detection algorithm 
will return maximum significance with a filter matched to the object size. 

Whether or not the galaxy is Gaussian, we can use equation 15411 to convert the logarithmic derivative expression, equation 
1511 1. into: 



dlnu 



2 ^ 
d, 



(56) 



since for near-circular Gaussian objects, e « e. (This follows because the smear polarizability P'^ is approximately the identity 
matrix in this case, as can be verified through direct computation using the formulas in Appendix|n|) The transfer hypertensor 
is then [see equation 1481 1: 



e ef}. 



(57) 



3.4 How Important is Shear Selection Bias? 

We can evaluate, at least at the order-of-magnitude level, the importance of the shear selection bias. Since ad is matched to 
the size of the object, the quantity 1 — a^/a^ will be near unity for well-resolved galaxies, decreasing to zero for point-like 
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objects. (That H = for pointlike objects is unsurprising because a point source remains a point source when sheared.) 
The T-factor is shown in Fig. |H1 The weighted ellipticities e are half of the adaptive ellipticities [equation ^] for near- 
circular Gaussian objects with a matched Gaussian weight function. Thus the e cir cular ellipticity has shape variance of 
approximately V ^ (0.32/2)^ 0.026, where 0.32 is the RMS ellipticity of galaxies ijMcKav et al.ll200lh . (This parameter 
will depend somewhat on the sample examined, so the value 0.32 should be thought of as a rough guide rather than a law of 
nature.) It thus appears that the shear selection bias will be of order several percent for well resolved galaxies, falling to zero 
in the limit of a poorly resolved object. 

A more quantitative approach to the problem is possible in the special case where the galaxy significance distribution is a 
power law and the shear estimator is constructed from an unweighted distribution of galaxies (i.e. w=constant). We suppose 
that the galaxy counts scale as TV oc The exponent A is determined by noting that v is a signal-to-noise ratio and hence is 
proportional to the galaxy's total fiux F divided by the square root of its post-PSF-convolved area A; the number-magnitude 
relation is roughly A'^ oc F~^. If the galaxies are poorly resolved, A is determined by the PSF and hence is constant, in which 
case A = 1. For well-resolved, nearby galaxies (redshift z ^ 1) we can assume a Euclidean universe, in which the usual inverse 
square laws F cc A cc apply; this gives F/VA oc 1/r cx iV"^/^, since the number of objects scales as . This gives 

A = 3. Thus we expect that A will probably lie in the range 1 < A < 3, with the distant galaxies of interest to cosmic shear 
measurements likely having A closer to 1. 

We next evaluate the average of equation (1571 . We define a shear selection bias resolution factor TZ = 1 — ag/a^, which is 
1 for well-resolved and for poorly-resolved objects. We also use that dN/dv x v^^^^ . The average transfer hypertensor is: 



, X^^-^-^P7^^(^;^o)e"e,d^ v-'-^Vnv;uo)Au 

Jo ^~^~^T^Av 

where TZeff is a weighted (by T) average of the shear selection bias resolution factor TZ over the detected galaxies: TZeff ~ 
{TZT}d/{T)d- Since T is a rapidly decreasing function of z/ (see Fig.|HJ, TZeff should be evaluated for galaxies near the detection 
threshold. 

At this point we note that the integrals in the numerator and denominator of equation 1581 diverge at the i? 
limit because 7^(0) > 0; see equation 1461 . However, this divergence must be unphysical since real source galaxy catalogs do 
not contain an infinity of faint objects. [The origin of the problem is that in deriving equation 14611 . we neglected the fact 
that objects overlap with each other.] We will thus remove this divergence by supposing that at very small i?, the detection 
probability V vanishes sufficiently rapidly for the integrals in equation l)58|l to converge; this yields (using integration by 
parts): 

r°^^-xav^- ■pJL(iy->')di, 

u~^~'^Vdu u^^^^Vdu 

Thus for cosmic shear measurements where galaxies are selected via a statistical significance cut, where V = 0.026 and A ~ 1, 
the shear selection bias causes an underestimate of the shear by 2.6TZe ff per cent. (It vanishes for the poorly resolved galaxies 
with Reff ~ 0.) This factor would be increased to roughly 8 per cent in the case of A = 3, Reff = 1, i.e. of low-redshift source 
galaxies for which the Euclidean source counts are more nearly valid and the near-threshold objects are well-resolved. Note 
that this result is independent of the significance threshold, at least when a power-law form for the galaxy count is assumed. 
Indeed, because it does not depend on the particular form for V, this result is still valid if only those source galaxies are 
selected that appear in multiple images of the same region of sky (this is the case in most current lensing surveys). It is also 
important to understand that shear selection bias is small not because the shear is small, but because the background galaxies 
are so nearly circular. If the galaxies are more elliptical in deeper surveys than in SDSS used here the bias will be larger. Note 
also that the mean RMS ellipticity, A and Reff are all likely to be redshift dependent and so will be the bias. 



3.5 Correcting for Shear Selection Bias 

While the above analysis serves to illustrate the order of magnitude of the shear selection bias problem, it may not be 
adequate to provide an analytical correction due to the large number of assumptions and approximations made. There are 
also contributions to the selection bias other than significance thresholds; for example, it is generally necessary to remove 
overlapping objects from the background galaxy catalog, but the overlap probability is a function of galaxy shape. There are 
four major approaches that could be taken to overcome shear selection bias: 

1. The bias could simply be ignored. 

2. The galaxy detection algorithm could be designed so as to eliminate, or at least substantially reduce, the bias. 
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3. The transfer hypertensor S could be computed and used as a correction to the shear polarizabihty or resolution factor. 

4. The bias could be determined from simulations. 

Of these, Method is appropriate when calibration errors of order several per cent are unimportant, for example in 
typical galaxy-galaxy lensing studies where statistical error dominates. However, it is clearly inapplicable to future cosmic 
shear surveys that aim for much higher accuracy. Method |5| requires a detection statistic that does not suffer from shear 
selection bias, i.e. one that is independent of the shear. Examples would be the galaxy's magnitude in a given band and the 
PSF-deconvoIved area. Unfortunately, these are not truly independent of the shear when noise is taken into account, and the 
"overlapping object" effect described above is not addressed. These effects are difficult to compute analytically, hence it is 
not clear whether we can verify that the residual errors in Method |5| are "small enough" without resorting to simulation. 

Method |3| requires that we compute Ep = \nV for each galaxy in our survey, then find H"/3 = e"i5/3, and compute the 
weighted average l^" i3)wd to determine the calibration error due to shear selection bias. The problem is thus one of finding 
Ep; unfortunately, Ep can only be determined from the intrinsic intensity / of the galaxy, not the observed image. Thus 
correction of cosmic shear measurements based on determination of H would require a careful analysis to ensure that noise 
and PSF effects do not invalidate the calculation. Method |4] the simulation, is straightforward albeit time-consuming. Its 
advantage is that it can accommodate arbitrarily complicated detection algorithms. The main weakness is that the underlying 
distribution of galaxies in the simulation may not be correct; however, the shear selection bias is only an effect of order several 
per cent anyway, so determination of the bias to many significant figures is unnecessary. 



4 CONCLUSION 

In this paper we have addressed several possible biases present in the data reduction of cosmic shear measurements. We focus 
on effects that affect the overall normalization only and cannot be detected using the usual methods of finding systematics, 
such as the tests on stars and E/B decomposition. We pay particular attention to the deviations from the elliptical Gaussian 
profile of the galaxy or the PSF, which is often the starting assumption of the data reduction methods. The results in Figs, 
mjishow that while for a Gaussian PSF most of the methods are successful in reconstructing the shear, the case of a more 
realistic non-Gaussian PSF is significantly more problematic except for the re-Gaussianization method introduced in this 
paper. This suggests that in existing surveys the shear calibration m ay b e inaccurate a t the 1 per cent level. Some methods, 
such as current implementation of BJ02 used in ljarvis et alJ l|200i) and lMcKav et alJ ll200lh . are particularly susceptibl e to 
this and can underestimate the l i near amplitude of fiuctuations by 10 per cent, while those bas ed on KSB95 (e.g. Bacon ct alj 
bond IVan Waerbeke et ahlEoonl . IPhodes et al.ll200ll. iHoekstra et allbooi iBrown et ai]l2002l) have a negative bias of up to 

5 per cent. It is interesting to note thar Tl£^i^etalT(|2002^ find a significantly (10-20 per cent) lower value for o"s than 
other surveys. Calibrations errors of the type discussed here can explain part of the difference, although other factors such 
as statistical errors, uncertainties in the redshift distribution of background galaxies, nonlinear evolution corrections of dark 
matter power spectrum and residual systematics (which are known to be present because their effect on B-type correlation is 
seen in all the surveys that have looked for it) also play a role. At the moment, residual systematics due to the PSF correction 
can be as great as 10 per cent in cts and are comparable to the statistical errors in current samples. 

There are other biases that one must worry about when accuracy on the order of one per cent is required. One such effect 
is shear selection bias which causes the galaxy selection to be correlated with the underlying shear and with the galaxy flux. 
This error is of order a few percent, but can reach up to roughly 8 per cent in some limits. Since the shear selection bias 
preferentially selects galaxies with shear orthogonal to cosmic shear the bias translates into an overall underestimate of the 
amplitude of fluctuations. Moreover, since for a flux limited survey the smaller and fainter galaxies are at a higher redshift 
the error in the amplitude may be correlated with redshift and can give an incorrect growth factor if not properly accounted 
for. Another effect that can create this is if the mean RMS ellipticity of the galaxies changes with redshift. 

Looking further ahead, weak lensing has the potential to break many of the degeneracies present in the cosmological 
parameters even after the cosmic microwave background has given us all of its information, especially if one can use the 
redshift information on background galaxies to extract the information on the growth factor It is an especially 

powerful probe of dark energy and massive neutrinos, since they both affect the growth factor at low redshifts probed by 
weak lensing. While the statistical precision of some of the proposed experiments such as LSST is impressive, the systematic 
effects must also be controlled at 1 per cent level. Current reduction methods are not sufficient for this purpose, but improved 
methods such as our re-Gaussianization method (which can be further generalized to include higher order terms) should be 
up to the task. Similarly, shear selection bias and other biases should be carefully simulated to estimate their importance as 
a function of galaxy detection significance, PSF etc. One notes that realistic simulations would allow us to determine both 
the PSF dilution calibration and the shea r select i on bia s , and indeed shear select ion bias is implicit ly taken into account 
in simulated observations such as those of lKaise3 ll2000l) . TErben et alJ ll200lh . and lBacon et al.1 ll200ll) . Future experiments 
will require that this line of work be extended, for example to include measured PSFs from the telescope and the correct 
distribution of galaxy shapes and fluxes as a function of redshift. 
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It is clear that the quahty of the data of the upcoming measurements will exceed the current generation of the reduction 
methods. However, this is not a cause for alarm, since none of the problems discussed here is fundamental in nature. We have 
proposed solutions for some of the problems, while others will require more detailed modeling using realistic simulations of the 
galaxies and their image processing through the atmosphere (if applicable), telescope, detector, and software. Given enough 
attention to these details weak lensing should fulfill its promise in the era of high precision cosmology. 
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APPENDIX A: CONVOLUTION WITH ELLIPTICAL LAGUERRE EXPANSION 

Here we examine the problem of convolving two functions / and g given in terms of their Laguerre expansions: 

f(^) = T,p,g apq1Ppq{x; Mf) 
g{x) = J2p,q bpqi'piix; Mg) 

The convolution is: 

h{x) ^ / f{x')g{x - x')cfx' = 2_]cpq''Ppqix;Mh); (A2) 

p,q 

we will consider here the case where M;, = M/ + Mg. This problem is easiest in Fourier space using the convolution theorem 
h{k) = 27r/(fc)g(fe); using equation 1101 and writing the result in terms of k = M^'^^fc yields: 

jJV -N'+N" 

pq p'q',p"q" 

We can simplify this by defining the matrices A'-*^' — MJ/'^M,^^'^^ and similarly A'^'. Then we use the result that any 2x2 
matrix, whether symmetric or not, can be written as a product of a rotation matrix, a diagonal matrix, and another rotation 
matrix (we will call this the RDR decomposition): 

^^/cos« -sinaWC U -sx "-^xA j^JC 0\ 

\^ sina cosa y \^ D y \^ sinx cosx / \ D J ^ ' 

where Re represents a rotation by angle 6, and C, D, and the rotation angles a and x can found explicitly from the 
rectangular-to-polar conversions: 

All + A22 ^{C + D) cos(a + x) 

A2i~ Ai2 = {C + D)8m{a + x) .... 
An - A22 = (C - 5) cos(q - x) ■ 
A21 + A12 = {C ~ D) sin(a - x) 

This decomposition is unique up to several discrete degeneracies: (1) interchange of C and D with phase shift of a by +7r/2 
and X by — 7r/2; (2) sign change of C and D with phase shift of a and x together by 7r/2; (3) phase shift of a and x together 
by tt; and (4) phase shift of a and x independently by multiples of 2-k. (In the special case that C — ±_D, there is an additional 
continuous degeneracy in that only a ± x, but not either angle individually, is determined.) If we decompose A'-^-* and A'"-*, 
then the condition Mj, = M/ -|- Mg yields that A'-'^'^^A'^' + a'^^'^A'"-' is the identity; this can be expressed in matrix form 
as: 



/ Cf cos^ X/ + Dj sin^ Xf (Dj - C/) sin X/ cos Xf \ 
\ iD) - C|) sin Xf cos Xf Cj sin^ Xf + D) cos^ x/ / 

/ C^cos^xs + ^sSin^X9 (^9 - f^s) sin Xg cos Xs \ ^ 1 0\ 
V (-^9 ~ ^l) sin X9 cos Xs Cl sin^ Xg + ^9 cos^ X9 / V ^ 1 / 

The trace of this equation yields + Cg + + — 2; the traceless-symmetric part yields: 
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cos2x/ \ , .^2 f^2^ ( cos 2X9 



which imphes that either (7^ = l)^ and Cg = D^, or that x/ and Xg differ by a multiple of 7r/2; we may, however, assume 
the latter case, since in the former case the aforementioned continuous degeneracy allows us to choose x/ = Xg- Indeed, 
the discrete degeneracy allows us to set Xf ~ Xg{= x)- From this we derive Cj + Cg = + Dg — 1. Finally, the positive 
determinant of A^-'^'' and A*-^' implies that CfDf and C'gDg are positive; the discrete degeneracy then allows us to choose Cf, 
Cg, Df, and Dg positive. 

Explicit expressions for a/, Og, x, Cf, Cg, Df, and Dg in terms of M/,g,h may be obtained by directly finding the square 
roots of these matrices, thereby determining A'-*^'^', and then substituting these matrices into equation IIA5fl . In order to 
construct the A's, it is necessary to determine the square root of a symmetric matrix. This is given by: 



j^l/2 ^ 1 f Mx:c + KUy 

a/m^^ + Myy + 2/D V ^-^yy + ^ 



(A8) 



M 



With this aside completed, we return to the determination of Cp, . Substituting the decomposition [equation IIA4H into 
equation 1A3II . and setting {U,V) — R_^k yields: 

■N -N' + N" , 

^^^==Cpge^™-A,,([/,y) = 2 ^p,|^,,p„,^„| V.'V."e'"""^e-"" -^A,,g,{CfU,DfV)A,ng„{CgU,DgV), (A9) 

pq p'q',p"q" 

where we have used the result that Apg(Refe) — e~""''Apg(fe). Then equation IIA9t is a linear equation for the Cpq, thus: 

Cp, = 20F.e-''"'^ ^ Z^;^'-^"'"(C/,5^)ap,,,V'9"e"""'"^e-™"'^«, (AlO) 

p'q' ,p"q" 

where the coefficients Zpg'' (Cf,Df) satisfy: 

•JV , ////__ jJV'+iV" 

^^^l,' ' {Cf,Df)kpg(U,V) = , „ „ „, ,,^ kp'q.[CfV,DfV)t.png„(CgU,bgV). (All) 



pq 

The Zpi^ ' {Cf,Df) can be calculated by noting that the right hand side of equation ijAllfl is a polynomial of order 
p' + q' + p" + q" in the two variables U and V. Since the Apq{U, V) with p + q < N form a basis for the (A'' + 1)(A^ + 2)/2 
dimensional vector space of polynomials in U and V of degree < A'', the sum on the left hand side of equation I|A11^ may 
be cut off at p + g < N. Equating coefficients of U^^V" on both sides then results in a finite system of linear equations 
for Zpg'' (Cf,Df), which is easily solved using standard linear algebra techniques. An alternative route is to exploit the 
orthonormality relations among the Laguerre basis functions to convert equation lAllll into: 

r.p'q',p"q",r i^^'+^^"-^^ KlgjU ,V) Kp. g, {C fU , D fV) Kp„ g„ [C gU , D gV)^- ^'"'^ dU AV 

^p\q\p'\q'\p"\q"\ j^, e-(u^+v')dUAV ' ^ ' 

this ratio of Gaussian integrals can be computed using Wick's theorem with the Feynman rules ([/) = iy) — (UV) = and 
iU'} = {V'} = 1. 

The Gaussian- integral method allows us to construct an explicit formula for Zpg'' ' (C/,-D/) in terms of finite sums. 
To do this, we expand the A: 

Apq{u,v) = ^ Xpg^^u^v", (A13) 

IJ,,v>0 

where: 

X - / Et(-l)''" %!(M-t)!(P-^-M+t)!t!(,-.-t)! s > 0, s an integer 
y otherwise 

and s = [p + q — n — p)/2. The sum in equation IIA14t is over all values of t such that s, jj. — t, t, and q — s — t axe all 
non-negative. The "selection rules" for X are as follows: Xp^^^ 7^ only if /i -|- is equal to p -f- q, or less than p -I- q by an 
even integer. Having constructed X, we evaluate the integral in equation IIA12t as: 

^'+^"'^ V- (2A-1)!!(2!>-1)!! 



^p\q\p'\q'lp"lq"\ ^ 2''+^^ 

J2 x;g^,Xp,g.^.^,Xp,.g..^..^„c^. Di {i^CfY" (i-DfY" . (ai5) 

IJ.-\- fi' -\- fj." —2 p, u-\-u' -{-h/" —2i> 
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We have computed several of the lower-order (up to degree 6) Z coefficients for the simulations in this paper; it is found 
that for these moments, it is practical to perform either the linear system solution [equation lAllH or the summation method 
[equation JA15H using symbolic manipulation software, and then write "hard-wired" C code to numerically evaluate the 
resulting expressions. Because the explicit polynomial form of the A functions [see equation ljA14^ contains many positive 
and negative terms for large p and q, it is likely that round-off errors would be an issue for our computational approach if we 
were to consider moments of sufficiently high order. It is unlikely that moments of very high order would be useful in weak 
lensing studies; nevertheless, if this turned out to be the case, or if the high-order Z coefficients were required for another 
application, it would be necessary to develop a better algorithm (likely involving a recursion relation) for computing them. 



APPENDIX B: LINEAR PSF CORRECTION 

We can derive a relation between e'^' and e'^', valid to first order in the departure of the galaxy and PSF from Gaussianity 
as follows. Suppose that the adaptive covariance of the intrinsic galaxy image is M/. Then the observed image I{x) is the 
convolution of / with the PSF g, and it can thus be expressed as a Laguerre expansion with covariance — ISAf + Mg. 
Using the convolution results of Appendix^ we find: 

I{x) = CV'p,(a;; M/) = Cp^i^p.ix; M,,), (Bl) 

pq pq 

where: 

Cpq=2^e---^ Y Z-;/^^''''\Cf,Df)e''"^'-fe-'^''-^blPX'l'^ (B2) 



and b^P , and b^^,} ,, are the Laguerre coefficients of / and g around M/ and Mg, respectively. (The rotation angles Xi Q^/i ^iid 
Ug and diagonal elements Of and Df can be calculated using the approximation M/ ~ M/ — Mg since corrections to this 
result will appear as second-order corrections to the final moments.) Since we are working to first-order in the non-Gaussianity 
of the galaxy and PSF, we may simplify equation IIB2t to give (for pq > 00): 



Coo 



.p'q'>00 p"i3">00 



(B3) 



where we have worked to first order in the /3's, and used the result that ^00'"° = 1. 

We next wish to relate the Cpq to the P^q^ . This can be done by use of the generators for point transformations provided 
by BJ02. To first order in the /3's and in Mi -Mh, and noting that = /322, we get [see BJ02 eqs. (6-45) and (6-46)]: 



(B4) 



Cll/Coo ~ fl 

C02/C00 ~ 
C20/C00 ^ 

Cpq/coo ~ I3iq^ pg/ 11, 02, 20 
where ^ is the dilation and rj is the shear required to transform from a covariance M/ to M?i. Explicitly: 

M. = (1 - 2,)Mi mr ( ll ) M^, (B5) 

where K and 3 denote real and imaginary parts, respectively. 

Assuming that we know the PSF, /S^??^,, can be determined. We further know that /9j(' = /Jgj' — ^20^ ~ 0- Now truncate 
equation ljB3|l with a finite set L of Laguerre modes, i.e. we sum only over p'q' G L. We then approximate the higher-order 
moments Cpq/coo ~ Ppq for pq £ L = L \ {00, 11,02,20} in accordance with equation IIB4ll . This gives us a system of linear 
equations that can be solved for the unknown /3pg''s: 

Y Zl/-''\C},D;)e-'"^'-ip\l,l,=e^"^^Pi[^ ^ Zf/' {C j , D p^^}^,, , (B6) 

p'q'GL p" q" £L 

for pq G L. The /3's so determined can be substituted into equation lB3t to yield ch/cqo and C02/C00, which then yield Mh 
via equations 1B41 and <B5ll . Finally, we can determine My = M^ — Mg and then determine the intrinsic ellipticity from 
equation @. 

The above procedure is somewhat involved, and we illustrate it with a simple example: that of L = {22}. That is, we will 
correct for the radial fourth moment or kurtosis of the PSF and intrinsic galaxy image. Obviously in ignoring the higher-order 
radial and angular moments we are making an approximation, however the 622 moment is present for non-Gaussian galaxy 
profiles such as the exponential or de Vaucouleurs profiles (and is positive in both cases). As we show in Section 12.51 the 
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higher-order moments can be important for some PSFs. The angular moments bpq for p ^ q are not present for pure radial 
profiles, i.e. profiles whose isophotes are concentric ellipses that share principal axes and axis ratios. If the angular moments, 
or the higher-order radial moments 633, 644, etc. are important, these would need to be included in L. 

Because of the shear and rotation invariance of the adaptive moments method, we can assume without loss of generality 



that the PSF has e'®' = 0, and that the measured object has e 



(I) 



and e 



{I) 



> 0, i.e. that it is elongated along the 



a:-axis. (Alternatively, if a rounding kernel has been applied to circularize the PSF, then we already have e'^' = 0, and no 
application of shear invariance is necessary.) We then have: 



It 





It 



and M/ 











Then to zeroeth order we determine: 



M 



/loth 



order = M/ - Mj = 



i[T,(l + eW)-r, 



(B7) 



(B8) 



___ 1/2 1/2 1/2 1/2 

Following the procedure outlined in Appendix|S| we compute the RDR decomposition of My M^^ and Mg M^^ ; this 
gives af — ag — X ~ 0, and: 



Co 
Da 



l-C,=T,/[T,(l + e(^')] 



Df=Tj[T,{l~e 



(B9) 



Next we must construct the system of equations, IIB6L which requires knowledge of the Z coefficients. The relevant Z 
coefficients for the intrinsic image are: 



5.22,00 _ 722,00 
^02 — ^20 



1 

2V2 



Zir^ = ±i^C^.+D]) + ^^C}D] 

Zir' = -\{C^- + D])-\C}D]+C] + D]- 



(BIO) 



and those for the PSF are: 



1 

2y/2 



700,22 



3,^ 



ZTi = -^K + Dt) - -ClDl + Ci, + Di,. (Bll) 
From these we construct the linear system, equation <B6ll : in this case it is only a 1 x 1 system and hence is easy to solve: 



-"22 ,722,00 



'-'22 



700,22 ^((9) 
^22 M22 



(B12) 



Then we can compute the correction M/ 1-^ due to non-Gaussianity using equations <B3t and (IB41 : 

— n In — 722,00 I 700,220(9) 

M — cii/coo - ^11 P22 +^11 P22 
r? = 2\/2co2/coo = 2\/2 



722,00 ^,(/) I 7' 
^02 P22 + ^' 



00,22 ^,(9) 
02 H22 



(BIS) 



(Note that in this case, r] is real. This is not necessarily the case if we use m 7^ modes in our basis L.) Finally, we compute 
M)h using equation IBSII and find that: 

1 / (l-2^i-r7)T,(l + eW)-r, 

2 I (\-2ii + ri)Ti(\-e^'^)-Tg 



Mo 



(B14) 



The ellipticity e'-'^ of the intrinsic galaxy image are easily computable from M/ using the covariance-ellipticity relations 
[equation Q]: 



M. 



if) 



(/) 



1 - 2^ - ri)Ti{l + e^") - T9] - [(1 ^2^I + ri)Ti{l - e^'>) - T,] 



[(1 - 2/x - ri)Ti[l + eW) - T,] + [(1 - 2p + ri)Ti[l - eW) - T,] 
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(l-2^-,7)(l + eW)-(l-2M + 77)(l-eW) _ -rj + {1 - 2^i)e^'^ ^^^^^ 



{l-2^I- J7)(l + eW) + (1 - 2m + J7)(l - eW) - 2Tg/Ti 1 - 2fi - rjeW - TJTi ' 
This is equivalent to BJ02's method, equation I14L except that the resolution factor is now: 

_ l-2M-77eW^r,/r, (,) 
-r, + (l-2M)ea) ' ■ 

If a shear has been applied to circularize the PSF, the inverse shear must now be applied to yield the physical e'-'^'. The 
transformation of ellipticities e e' under a shear {5+,8^) (i.e. a shear that maps a circle to an ellipse of ellipticity 5) is 
given by [see, e.g. BJ02 eq. (2-13)]: 



e 



1 

l + 5^e^ + d"x ex 



+ 5+ + %(1 - \/r^)(ex<5+ - e+5x) 



+ 5x + ^(1 - VT~S^){e+5^ - ey,5+) 



(BIT) 



where 5'^ — 5\ + 5%^; the transformation of traces can be determined from conservation of area, 

[Warning: the "shear composition" operator e' = 5 © e as defined by BJ02 and in equation IB 1711 does not form a group 

because it is not associative. This problem arises because the shear e followed by the shear 5 does not have the same effect 

as the shear e' (they differ by an initial rotation). For this reason, it is probably better to think of © as representing the 

transformation law for ellipticities under shear, rather than as a composition operator.] 

We can thus list the steps involved in PSF correction using the moments L = {22}. If the PSF has been circularized 

using a rounding kernel, the steps marked with a * can be skipped. 

• Compute the adaptive covariances M/ and Mg for the galaxy and PSF, and the corresponding radial fourth moments (3':^2 
and4^\ 

• ★ Apply a shear to circularize the PSF, i.e. use equation HBlTt with 5 — — e'®'. The radial fourth moments are unchanged 
by this shear. 

• Apply a rotation so that the measured galaxy has e!^' > and e'^^' = 0. 

• Compute the dilation /i and shear rj in accordance with equations lB7t through JB13II . 

• Compute the intrinsic ellipticity of the galaxy in the sheared rotated coordinate system using equation (IBlSt . 

• Undo the rotation. 

• ★ Shear to return the ellipticity to the unsheared unrotated coordinate system. The shear is undone using equation lB17t 
i.e. the exact opposite shear that was applied initially. 



APPENDIX C: THE KSB95 METHOD IN LAGUERRE COEFFICIENTS 

Here we describe our im plementation of the KSB95 shear estimator, updated to include the pre-seeing shear polarizability 
iLuppino fc Kaise3ll997^ . For simplicity, we will assume a Gaussian weight function. We begin by computing the Laguerre 
expansion [see equation llH l of the measured object I and PSF g to second order with a circular covariance matrix: 

bl,[^ = —L=f Iix)e-^''"+y"'>/^"\;,{x/a,y/a)d^x, (CI) 
a^^n-p]q] J ^2 

where a is the width of the weight function. Then KSB's ellipticity (which in general is distinct from the adaptive 
ellipticity!) is given by: 

K , ■ K \/2620 m u , u /r^n\ 

e+ + lex = — — , Ttu = boo + bn. (C2) 

J- W 

Then the shear polarizabihty hypertensor is given by: 

[0./3 = [X*],^-[e^^%[e>, (C3) 



where 

^shi boo ~ b22 c , %/6 / 5R&40 9^40 \ /^.N 



and 



e+ +iex = . (C5) 

J- w 

The KSB/LK method also requires the smear polarizability hypertensor, given by: 

[P='^./3 = [X^n^fi - [e^'Ue'^h, (C6) 
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where 

rv-smi 600 — 26ll + 622 e , ( 5R&40 Q'&40 1 

]"^ = ^"''+47^U&40 ' ^^^^ 



and 

V^&20 — V6&31 

4r„ 



(C8) 



In these equations, the X contributions to the polarizabihties are due to the effect of shear on the quadrupole moments &20 
whereas the e^e^''''^™ contributions are due to the effect of the shear on the weighted trace Tw When the X hypertensors 
are expressed in terms of the Laguerre moments, the decomposition into spin-0 (scalar) and spin-4 (traceless-symmetric) 
com ponents is ma nifest . 

iLuppino fc Kaiser 1 199^ then defines a pre-seeing shear polarizability hypertensor to take account of the finite PSF size: 

[P'U = [P'''{I)U - (C9) 

Our only remaining issue is the choice of the width a of the weight function. We use a weight function scale set by the 
requirement &11 = 0, i.e. we use an adaptive circular variance. This is appropriate, at least for the simple case of a circular 
Gaussian object, for maximum significance detection with a Gaussian filter of variable width (see the discussion in Section 
ESI- FoUowine- Hockstra et al.. L199&) . we use the same a for the PSF as for the object. 
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